THE SPECTRAL BASIS AND RATIONAL INTERPOLATION* 



GARRET SOBCZYK t 

Abstract. The Euclidean Algorithm is the often forgotten key to rational approximation tech- 
niques, including Taylor, Lagrange, Hermite, osculating, cubic spline, Chebyshev, Fade and other 
interpolation schemes. A unified view of these various interpolation techniques is eloquently ex- 
pressed in terms of the concept of the spectral basis of a factor ring of polynomials. When these 
methods are applied to the minimal polynomial of a matrix, they give a family of rational forms of 
functions of that matrix. 

Key words. Euclidean algorithm, cubic spline, interpolation, rational interpolation, spectral 
basis. 

AMS subject classifications. 13F10, 13F20, 15A24, 41A10, 41A15, 41A20, 41A21, 65D05, 
65D07, 65D17 

1. The Euclidean Algorithm and Spectral Basis. The euclidean algorithm 
has many important well-known consequences in number theory, algebra and anal- 
ysis. In the spirit of pOj, we are mainly interested here in some of its consequences 
regarding interpolation of functions over the real or complex numbers. Let lR[x] and 
(Z|z] denote the rings of real- valued and complex- valued polynomials over the real and 
complex number fields IR and W, respectively. Whereas we state our results in terms 
of polynomials over the field of real numbers M, all of the results are equally valid for 
polynomials over (U. 

Let h(x) denote the monic polynomial defined by 



(1.1) h{x) = '[[{x-x,y 



where {xi, . . . , x,.} are the distinct real roots of h{x) with multiplicities {mi, . . . , m^}, 
respectively. Let f{x) E lR[x] be any polynomial in M[x]. Then the euclidean algo- 
rithm simply tells us that there will always exist polynomials g{x) and a remainder 
r{x) such that 

(1.2) f{x)^g{x)hix)+r{x) 

where < deg(r(x)) < deg{h{x)) or r{x) = 0. When equation p.2|l holds, we say that 

f{x) = r(x) mod(/i) where h = h{x), or more concisely, that f{x) = r{x). We denote 
the ring of all real polynomials modulo h{x) by IR[x\h. In the terminology of factor 
rings, Wt[x\h=R[x\/ < h{x) >, meaning that lR[x]h is isomorphic to the factor ring 
M[x]/ < h{x) > of M[x] generated by the principal ideal < h{x) >, O p. 266]. Thus, 
iR[a;]h has the structure of a ring with addition and multiplication of polynomials 
defined modulo h(x). 

By the standard basis of residue classes of ]R[x]h we mean 

(1.3) Bh = {l,x,x^,...,x"'-^}, 
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where m = deg{h) — mi + . . . + rrij.. However, calculations in iR[a;]/, are much more 
simply carried out by appealing to the special properties of the spectral basis [7], 
[S]. The spectral basis Sh of ]R[x]h consists of idempotents Si — Si{x), and nilpotents 
Qi = qi{x) and their powers = qf{x), 

/I /i\ c t mi— 1 mo — 1 m^ — It 

(1.4) Sh^{si,qi,...,qi^ ,S2,g2,...,g2 , ■ ■ ■ , Sr, qr: ■ ■ ■ , qr }, 

which satisfy the following properties under addition and multiplication in lR[x]h 
modulo h{x): 

Property 1. si + S2 + ■ ■ ■ + Sr = I, and SiSj = SijSi for i, j — 1, . . . ,r where Sij — 
for i =/= j and Stj = 1 for i = j. 

Property 2. qiSi = qi, and g™'^^ ^ mod{h) but g™' = 0, for i = 1, . . . , r. 

Property 3. For each /(x) e iR[a;], /(a;)si = (/(x)mod[(a;— a;i)™'])si for i = 1, . . . , r. 

Property 1, shows that the Si{x) are mutually annihilating idempotents which 
partition unity. Property 2, shows that Si acts as an identity element when multiplied 
by the nilpotent qi{x), and that qi(x) is a nilpotent of index in lR[x\h- Property 
3, shows that for each polynomial f{x) € IR[x\, Si{x) acts as the projection of /(x) 
onto the ring of polynomials -/R[x](j._2,.)»ni modulo (x — x^)™'. 

From these three algebraic properties, we can explicitly solve for the polynomials 
that make up the spectral basis. For each i = 1, . . . , r, define hi = hi{x) — h{x)/{x — 
Xi)"^' . Using Properties 1 and 3, we find that 

(1.5) hiSi = hi, 

and since hi(xi) ^ 0, it follows that /i,^^ 7^ mod (x — Xi)™'. Multiplying both sides 
of equation (|1.5() by /i^^niod (x — xj™', we find that 

(1.6) s,(x) = (/i-^mod (x - x,r')h,{x), 

which gives an explicit solution by taking the first terms of the Taylor series for 
h~^ around the point x = Xi. Having found the idempotents s^, the corresponding 
nilpotents g^, and their powers, are specified by 

q^ : = (x - Xi)^Sr = (x - x^f {h^r^mod (x - Xi)"''^)h^{x) 

(1.7) = (/i- Vod (x - x,)'"'-'=)/i,(x), 

for k = 0, 1, . . . , mi — 1. Note that for fc = 0, we get the correct convention that 

h 

The transition from the standard basis H1.3fl to the spectral basis (|1.4|l is accom- 
plished by first noting that 

r r r 

(1.8) X = ^ XSi = E((^ ~ ^ ^^^^^ ^ E^^* ^ 

i—l i—1 i—1 

from which it follows that 

r r nii — 1 / , \ 

1=1 i=l j=0 ^ ^ 

for fc = 0, 1, . . . , m — 1, as easily follows by referring to the properties of the spectral 
basis [H]. 
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2. Rational Interpolation. Let f{x) be a real-valued lunction which is contin- 
uous and has derivatives to the orders {nii — 1, ... , rrij. — 1} at the respective points 
{xi, . . . , Xr}, where {xi, . . . , Xr} are the distinct real roots ol h{x) with multiplicities 
{mi, . . . , rrir} as was defined in (|l.l|l of the previous section. 

The function f{x) of the real variable x can be extended to a function of the 
variable x e lR[x]h by simply substituting p. 8(1 into /(x), getting 

r r 

(2.1) f{x) := /(^(x, + qi)Si) ^ ^ /(x, g.)^,. 

i=l i=l 

If we now expand f{xi + g^) in a Taylor series about x = Xi^ we get the desired 
expression 

mi — 1 . 

(2.2) /(x. + gO ^ 5] 



fc=0 



where as usual, /'■'^•'(xi) = ^^/(a;)U=xi- Although we have derived equations H2.1|l 
and (|2.2|l modulo h{x) from the basic properties of the spectral basis (|1.4|l . we could 
equally well have taken 1(2. l|l and (|2.2|l to be the definition of /(a;) modulo /i(a;). 
The interpolation polynomial 

r m^— 1 

(2.3) g{x) = fix) mod M^) = Et E fci/^'^^O^f 

i=l fc=o 

is called the Birkhoff or osculating interpolation polynomial of f{x) with respect to 
h{x). In the special case that mi = ••• = nir = 1, g{x) is called the Lagrange 
interpolation polynomial of f{x), and in the special case when mi = . . . = rrir = 2, 
g{x) is called the Hermite interpolation polynomial of /(x), 5, pps.278,287], PI p. 52]. 
When r = 1, g{x) reduces to the first mi — 1 terms of the Taylor series of f{x) about 

X = Xi. 

Rational interpolation also takes an equally eloquent form when expressed in 
terms of the spectral basis p.58],|2. Let a{x) = X^i'io^ a-ix"^, and b{x) = X^Jl^^ ^j^^ 
be polynomials over the real numbers IR. We say that 

2.4 g{x) = 

b{x) 

is a rational interpolate of f{x) at the points (nodes) {xi, . . . , Xr} with multiplicities 
{mi, . . .,mr} if 

(2.5) f{x)b{x)-a{x) ^0. 

The usual way of defining rational interpolation involves the solution of a system of 
linear equations for coefficients of the polynomials a{x) and b{x). Our definition is 
simpler and more direct in that it only requires that the modular relation (|2.5|) holds. 
When b{x) has no common zeros with h{x), the equation (|2.5|) is equivalent to 

(2.6) fix) ^ a{x)b{x)-\ 

Essentially, each choice of the shape polynomial b{x) in (|2.5|) and H2.6|l determines 
a different rational interpolation g{x) of f{x), Because of the homogeneous 
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nature of the rational interpolate (|2.4|l . we can require that b{xi,) = 1 at some point 
Xf. which is chosen not to be one of the roots Xi of h{x) nor a zero of f{x). Whereas 
equation (|2.5|l is defined niodulo(/i(x)) using (|2.1|l and H2.2|l . equation H2.4|l is an 
ordinary equality. 

Chebyshev and other kinds of rational interpolation are defined simply by re- 
placing the powers of x'' of x in a{x) and b{x) by the corresponding Chebyshev or 
other sets of orthogonal polynomials of the same order. A comprehensive study of 
the algebraic structure of rational functions has been undertaken by Luis Verde-Star 
in a series of papers ^1 1111 [T^ . Examples of the various kinds of interpolation will 
be given in Section 3. 

Cubic spline interpolation can be defined parametrically in terms of the spectral 
basis ^2,2 of lR[t]h for h = h{t) = t^{t - 1)^. Using the formulas ((T^ and ((T~7|) from 
the previous section, we find that 

(2.7) 52,2 = {si = {2t + \){t - 1)2, 91 = t{t ~ If ,32 = (3 - 2t)e,q2 = f{t - 1)}. 

The piecewise natural cubic spline {51 (ii), 52(^2), ■ • ■ j <?fc-i(^fc-i)}, for < < 1 and 
A; > 3, connecting the successive points {xi, X2, . . . , x^} in 5?", is defined by 

(2.8) gi{ti) = XjSi(t^) -f ViqiiU) + x,j+iS2(ti) + Vi+iq2{U), 
with the requirements that 

(2.9) 5^(0) = = 5^1 (1) and 5^(1) = 5-'+i(0) 

for i — 1, . . . , fc— 2. Taking the second derivatives of gi(ti), and evaluating at = 0, 1 
gives 

(2.10) 5^(0) = 6(x,+i - x,) - 4v, - 2v,+i, 
and 

(2.11) <?f (1) = ~6(x,+i - X,:) + 2v, + 4v,+i. 

The resulting k linear vector equations are uniquely solved for the fc-unknown tangent 
vectors Vi,V2, . . . ,Vfc. 

If, instead of the requirement H2.9|l . the vectors vi and v^^ are taken as given, and 
the remaining (fc — 2) linear vector equations 

(2.12) 5-'(l)=.9:'+i(0), 

for i = 1, . . . , fc — 2, are uniquely solved for the (fc — 2) unknown tangent vectors 
V2, V2, . . . , Vfc_i, the resulting solution is called the bounded cubic spline, p. 93]. 

Various kinds of rational cubic splines can also be easily constructed by replacing 
the spectral basis 52,2 in (|2.8|) by a rational spectral basis of the form 

(2.13) 7?.2,2(6) = {Srl, grl, Sr2, qr2} 

for 

, , b{t)si{t) mod(/i(t)) b{t)qi{t) mod{h{t)) 

'^'^'^ ^ m ' ^ ^ W) 
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and 

b{t)s2{t) mod{h{t)) b{t)q2{t) mod{h{t)) 

m ' ''■'^ ^ ^ m ' 

where b = b{t) = l + bit + b2t^ + bQt^ and 6(1) 7^ 0. Clearly, the rational spectral basis 
TZ2,2{b) reduces to the ordinary spectral basis ^2,2, given in H2.7|l . for b — b{t) = 1. Of 
course, when using the rational spectral basis, the second derivatives g'/{0) and g'/{l), 
given in l|2.1()|l and H2.11(l . must be recalculated. 

3. Circles. A circle and other conies are good geometric figures on which to 
carry out interpolation experiments. We derive here several approximations for the 
unit circle, centered at the origin, using rational spectral bases. 

The rational spectral basis Sh = {sri, Sr2, Srs} for 6=1 + bit + 62^'^ and h{t) = 
{t + l)t{t - 1), is defined by 

_ \{l-bi+b2)t{t-l) (t-l)(f+l) 
l + bit + b2t^ l + 6it + 62t2' 

and 

_ ^(l + 6i+62)t(t + l) 
"^"^ " 1 + bit + 62i2 

We wish to optimize (in the sense of least squares) the choice of 61 and 62 so that 
the interpolating curve g{t) = /(— l)sri + /(0)sr2 + /(l)sr3 to the semi-circle f{t) = 
( cos(7ri/2), sin(7rt/2) ), for —1 < t < 1 is as good as possible. The values of 61 and 
62 can easily be found by requiring that 5(1/2) • 5(1/2) — 1 ~ g{—l/2) ■ g(— 1/2), 
giving the values 61 = and 62 = ±1. The value 62 = —1, gives the single point (1, 0), 
whereas 62 = 1 gives the well-known parameterization of the circle g{t) = {j^, TTT^)' 
A family of rational approximations to the quarter unit circle through the nodes 
(1,0) and (0,1), with the initial and terminal tangent velocity vectors (0, t;) and 
(—11,0), is specified by 

g{t) = (1, 0)S^1 + (0, V)qrl + (0, 1)S,.2 -f {-V, 0)qr2, 

in the rational spectral basis H2.13|l . One popular construction of the circle is based 
on NURBS (nonuniform rational B-splines) p. 110]. Letting 61 = —2 -I- ^2 = —62, 
the choice v — \f2 precisely eliminates the t^ term in the numerator, and gives the 
nurb parameterization 

git) = (1, 0)s,i + (0, V2)q,i + (0, l)s.2 + (-^2, 0)(j.2 



1 + (-2 + V2)t + (1 - \/2)i2 V2t+(1-V2)t2 



1 + (-2 + V2)t + (2 - \/2)i2 ' 1 + (-2 + V2)t + (2 - \/2t2 

A quite different perfect parameterization of the unit circle through the inter- 
polation points (1,0) and (0,1), and taking the initial tangent vector at (1,0) to be 
(0,7r/2), can be derived by using the rational spectral basis of 

52,1 = {si = -{t + l)it- l),qi = + t)t, S2 = t^} 
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Fig. 3.2. Both sine and cosine curves are shown together with their approximations. 

for h{t) ^t^{t-l). We find that 

(1, 0)si + (&i(l, 0) + (0, f )) qi + (0, 1)(1 + h+ 62)S2 
5W = ^ 1 ^, ,i, ,2 

l+bit + 02* 

Optimizing bi and 62, we find that 

8 + 4(7r - 4)t - 4(7r - 2)*^ 47rt + (tt - A)TTt'^ 



9{t) 



+ 4(7r - 4)i + (7r2 - 47r + 8)*^ ' 8 + 4(7r - 4)i + (tt^ - 47r + S)*^ 



Whereas the above parametcrizations give perfect circles, there are many other 
parameterizations that are interesting. For example, consider the family of approxi- 
mations to the miit semicircle through the points (1, 0) and (—1, 0), given by 

git) = (1, Q)sri + (0, v)qri + (-1, 0)s,p2 + (0, -v)qri 

employing the rational spectral basis of the kind (|2.13|l . Choosing h = 1 — t + t^ , and 
V = i gives a very good approximation to the unit semicircle for < i < 1, 

, /I -t - 3f + 2t^ ~3(t-l)t\^, 
(3.1) 9{t)^[ l-t + t^ ' l-t + t^ j=(cos7rt,sm7r0. 

with a least square error less than .000071, see figure 3.1 It is interesting to note 
that this parameterization also gives a good approximation to cos irt and sin irt for 
< t < 1, see figure 3.2 The series expansions for the approximations to cosTri and 
siuTrt are 

oc 

cosTTtSl - - 21^ + Y,{-lf+^[2t'''+^ + 4t3fe+2 + 2t"'+^] 

k=l 

and 

simrt^St + 3^(-l)'^[t3'= + t^''+^], 

k=l 

which are interesting in their own right. 
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4. Matrices. Let A be any n x n matrix over a field JC. The field JC may be 
the real or complex numbers, or even a finite Galois field. It is well known that every 
matrix satisfies it's characteristic polynomial, defined by 

ip{x) = det(2;/ - A), 

where / is the identity n x n matrix [HI 21 ■ The theory of a spectral basis is directly 
applicable to a matrix whenever the characteristic polynomial is of the form ip{x) = 
nl=i(^ — li)"' for distinct roots Xi G /C. When applying the spectral basis to a 
matrix A, it is better to use the minimal polynomial ipi^) = 111=1 ~ 2;^)™% where 
1 < mi < Ui for each i — l,...,r. The minimal polynomial of the matrix A is 
defined by the condition that it is unique monic polynomial of least degree for which 
^{A) = 0. 

Any of the interpolation formulas, developed in terms of the spectral basis in the 
previous sections, apply immediately to the matrix A, provided that h{x) = ^{x). 
This is because the relationship that h{x) = mod h is precisely reflected in the 
condition that ^p{A) = for the minimal polynomial tp of the matrix A. Thus, the 
spectral form of the matrix A is given by simply replacing x by the matrix A in (|1.8|l . 
getting 

r r 

(4.1) A^"^XiS,{A) + qi{A) ^"^XiS^ + Qi, 

1=1 i=l 

for Si = Si{A) and Qi = qi{A). The matrices {Si,Qi} of the spectral basis satisfy ex- 
actly the same rules under matrix addition and multiplication as does the polynomials 
{si,qi} of the spectral basis modulo h{x), [7|. 
For example, the matrix 



A- 



-2 -1 0^ 

-15 11 

-6 -5 -1 

-8 -10 -3 0/ 



has both characteristic and minimal polynomials /i = t^(i — 1). As a consequence, we 
can apply all the above interpolation formulas, found for the rational spectral basis 
of h{x), without modification. Using (|2.3|) . 12. 7|) . H3.1|l and H4.1|) . we find that 



cos ttA = / - + 4A^ 



I-A~3A^ 
I-A + A^ 



2A^ 



I 



1 2 

-10 -3 



10 
V 20 



8 

14 



-10 


6 


26 


-14 


-45 


26 


-82 


47 



and 



sinTrA = A - A^ ^ 



2 -3^2 ZA 



I-A + A^ 




As a check, we calculate cos^ tvA -\- siv? nA = I as expected. 
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